#ifndef ALGORITHM_H
#define ALGORITHM_H

/*
  算法模块

  作者：李江军
  时间：2020-11-20
  说明：算法模块开始使用于数据采集软件项目中，该软件隶属北京旋极信息科技股份有限公司

*/

#include "algorithm_global.h"
//#include "../include/commdefine.h"
#include <complex>
#include <QObject>
#include <QVector>
#include <QtMath>



/*
 * 频谱分析函数,快速傅里叶变换
 *
 * 作者：李江军
 * 时间：2021-07-28
 * 功能：定义滤波使用参数
 */
// 滤波器类型
enum FilterType
{
    LOWPASS,        // 低通
    HIGHPASS,       // 高通
    BANDPASS,       // 带通
    BANDSTOP        // 带阻
};
// 窗口类型
enum FilterWindowType
{
    RECTANGLE,      // 矩形
    TUKEY,          // 图基
    TRIANGLE,       // 三角
    HANN,           // 汉宁
    HANNING,        // 汉明
    BRACKMAN,       // 布拉克曼
    KAISER          // 凯塞
};






/*
 * 频谱分析函数,快速傅里叶变换
 *
 * 作者：李江军
 * 时间：2020-11-23
 * 功能：输入原始波形，输出频谱波形
 * 参数：in：输入参数，输入原始波形，输入数据个数必须为 2^n 的整数倍，否则无法计算
 *      out：输出参数，输出频谱波形，输出长度与输入长度相同
 */
template <class _Ty>
void FourierTransform(const QVector<std::complex<_Ty> >& input, QVector<std::complex<_Ty> >& output)
{
    int n = input.size();		//输入的点数，n>0
    int it, k(0);

    for(it = 1; ; it *= 2)	//求k值, 要满足n等于2的k次方
    {
        if(it == n) break;
        k++;
        if(k > 64) return;
    }

    output = input;
    QVector<std::complex<_Ty> > temp(n);
    for(it = 0; it < n; it++) {
        int m = it, is(0);
        for(int i = 0; i < k; i++) {
            int j = m / 2;
            is = 2 * is + (m - 2 * j);
            m = j;
        }
        temp[it] = output[is];
    }

    output[0] = std::complex<_Ty>(1.0, 0.0);

    _Ty p = 2 * M_PI / (1.0 * n);

    output[1] = std::complex<_Ty>(cos(p), -sin(p));

    for(int i = 2; i < n; i ++) {
        p = output[i-1].real() * output[1].real();

        _Ty q = output[i-1].imag() * output[1].imag();
        _Ty s = (output[i-1].real() + output[i-1].imag()) * (output[1].real() + output[1].imag());

        output[i] = std::complex<_Ty>(p - q, s - p - q);
    }
    for(it = 0; it<n-1; it = it + 2) {
        _Ty vr=temp[it].real();
        _Ty vi=temp[it].imag();

        temp[it] = std::complex<_Ty>(vr + temp[it+1].real(), vi + temp[it+1].imag());
        temp[it+1] = std::complex<_Ty>(vr - temp[it+1].real(), vi - temp[it+1].imag());
    }

    int m=n/2;
    int nv(2);

    for(int l0 = k-2; l0>=0; l0--) {
        m = m / 2;
        nv = 2 * nv;
        for(it = 0; it <= (m - 1) * nv; it = it + nv)
          for(int j = 0; j < (nv / 2); j ++)
          {
              p = output[m*j].real() * temp[it+j+nv/2].real();
              _Ty q = output[m*j].imag() * temp[it+j+nv/2].imag();
              _Ty s = output[m*j].real() + output[m*j].imag();

              s *= (temp[it+j+nv/2].real() + temp[it+j+nv/2].imag());
              _Ty poddr = p - q;
              _Ty poddi = s - p - q;

              temp[it+j+nv/2] = std::complex<_Ty>(temp[it+j].real() - poddr, temp[it+j].imag() - poddi);
              temp[it+j] = std::complex<_Ty>(temp[it+j].real() + poddr, temp[it+j].imag() + poddi);
          }
    }

    for(int i=0; i < n; i++)
    {
        output[i] = std::complex<_Ty>(sqrt(temp[i].real() * temp[i].real() + temp[i].imag() * temp[i].imag()),
                                   output[i].imag());
        if(qAbs(temp[i].real()) < 0.000001 * qAbs(temp[i].imag()))
        {
            if((temp[i].imag() * temp[i].real()) > 0)
                output[i] = std::complex<_Ty>(output[i].real(), 90.0);
            else
                output[i] = std::complex<_Ty>(output[i].real(), -90.0);
        }
        else
            output[i] = std::complex<_Ty>(output[i].real(), (atan(temp[i].imag() / temp[i].real()) * 360.0 / M_PI * 2));
    }
}

/*
 * 频谱函数
 *
 * 作者：李江军
 * 时间：2021-08-19
 * 功能：输入原始波形，输出频谱波形
 * 参数：src：输入参数，输入原始波形
 *      values：输出参数，输出频谱波形
 */
extern void ALGORITHMSHARED_EXPORT frequencySpectrum(const QVector<double> &src, QVector<double> &values);





class ALGORITHMSHARED_EXPORT Algorithm : public QObject
{
    Q_OBJECT
public:
    Algorithm(QObject *parent = Q_NULLPTR);

    /*
     * 频谱分析函数
     *
     * 作者：李江军
     * 时间：2020-11-20
     * 功能：输入原始波形，输出频谱波形
     * 参数：in：输入参数，输入原始波形
     *      out：输出参数，输出频谱波形
     */
    void FFT(const QVector<std::complex<double>> &in, QVector<std::complex<double>> &out);

    /*
     * 作者：李江军
     * 时间：2020-11-20
     * 功能：获取一个长度的最接近 2^n 的值
     * 参数：in：实际长度
     * 返回值：接近 2^n 的值
     */
    int powerValue(int len);

    int get_computation_layers(int num);

};


class ALGORITHMSHARED_EXPORT Filter
{
public:
    /* fType:过滤器类型
     * wType:窗口类型
     * scala:阶数
     * freq:频率
     * lfreq:下限频率
     * hfreq:上限频率
     */
    Filter(FilterType fType, FilterWindowType wType, int winLen, qreal freq, qreal lfreq, qreal hfreq);
    ~Filter();

    void setFilterType(FilterType fType);
    void setFilterWinType(FilterWindowType wType);
    void setWinLength(int len);
    void setFreqency(int freq);
    void setFreqUpper(int fUpper);
    void setFreqLower(int fLower);
    void setFilterParams(FilterType fType, FilterWindowType wType, int scala, qreal freq, qreal lfreq, qreal hfreq);
    // 获得滤波器系数后,对输入信号和系数做卷积积分,得到滤波器输出,d:原始数据
    void Convolution(const QVector<qreal> &d);

protected:
    void calcImpactData();
    qreal calcByWindowType(int WinType, int Scala, int Index, qreal Beta);
    qreal calcKaisar(int Scala, int Index, qreal Beta);
    qreal bessel(qreal Beta);
public:
    QVector<qreal>      ImpactData;     // 冲击响应输出信号
    QVector<qreal>      OutData;        // 滤波后输出数据

    qreal               lower;
    qreal               upper;
private:
    FilterType          Type;
    FilterWindowType    WinType;
    int                 Scala;
    qreal               Freq;
    qreal               lFreq;
    qreal               hFreq;
    QVector<qreal>      RawData;

    /**
      用于连续滤波，完成一组滤波后保存最后一组数据
     */
    QList<qreal>      m_preFilter;
};

#endif // ALGORITHM_H
